Structure of the Galactic Stellar Halo Prior to Disk Formation 



Masashi Chiba 

National Astronomical Observatory, Mitaka, Tokyo 181-8588, Japan 
email : cliibams@gala.mtk .nao.ac.jp 

Timothy C. Beers 

■ Department of Physics & Astronomy, Michigan State University, E. Lansing, MI 48824 

■ email: beers@pa.msu.edu 
O " 
(N . 

. . ABSTRACT 
> 

o 

I We develop a method for recovering the global density distribution of the ancient 

\^ ' Galactic stellar halo prior to disk formation, based on the present orbits of metal-poor 

stars observed in the solar neighborhood. The method relies on the adiabatic invariance 
^ ■ of the action integrals of motion for the halo population during the slow accumulation 

. of a disk component, subsequent to earlier halo formation. The method is then applied 

to a sample of local stars with [Fe/H]< —1.5, likely to be dominated by the halo 
component, taken from Beers et al.'s recently revised and supplemented catalog of 
I metal-poor stars selected without kinematic bias. We find that even if the Galactic 

' potential is made spherical by removing the disk component in an adiabatic manner, 

the halo density distribution in the inner halo region (i? < 15 kpc) remains moderately 
flattened, with axial ratio of about 0.8 for stars in the abundance range [Fe/H]< —1.8 
and about 0.7 for the more metal-rich interval —1.8 <[Fe/H]< —1.5. The outer 



o 



, halo remains spherical for both abundance intervals. We also find that this initial 

^ ' flattening of the inner halo is caused by the anisotropic velocity dispersions of the 

halo stars. These results suggest that the two-component nature of the present-day 
^ ' stellar halo, characterized by a highly flattened inner halo and nearly spherical outer 

■ halo, is a consequence of both an initially two-component density distribution of the 

halo (perhaps a signature of dissipative halo formation) and of the adiabatic flattening 
of the inner part by later disk formation. Further implications of our results for 
the formation of the Galaxy are also discussed, in particular in the context of the 
hierarchical clustering scenario of galaxy formation. 

Subject headings: Galaxy: evolution - Galaxy: halo — Galaxy: abundances — Stars: 
Population II 



1. Introduction 



The currently favored cold dark matter (CDM) theory of galaxy formation postulates that 
the formation of a massive spiral galaxy like our own is a consequence of the hierarchical assembly 
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of subgalactic dark halos, and the subsequent accretion of cooled baryonic gas in a virialized, 
galaxy-scale dark halo (e.g., Peacock 1999). Numerical studies based on this picture are able to, 
at least qualitatively, reproduce the characteristic features of a disk galaxy - the massive dark 
halo, the stellar halo, and the stellar disk components (e.g., Steinmetz &; Miiller 1995; Bekki & 
Chiba 2000; Navarro & Steinmetz 2000), though difficulties are still encountered in the details. 
For example, the simulations conducted to date do not adequately account for the size of the 
disk component and the number of satellite galaxies (Navarro, Prenk, k, White 1995; Moore et al. 
1999; Klypin et al. 1999). 

The CDM hierarchical model may be regarded, in its essence, as a modern generalization of 
the classical Searle & Zinn (1978) hypothesis for the formation of the Galactic stellar halo. To 
explain a large inferred spread in the ages of globular clusters, and the lack of a spatial gradient 
in their metal abundances, Searle &: Zinn argued that the halo component may have experienced 
prolonged, chaotic accretion of subgalactic fragments, as opposed to the rapid, monolithic collapse 
proposed by Eggen, Lynden-Bell, & Sandage (1962). Recent discoveries of halo substructures in 
Galactic phase space (Majewski, Munn, & Hawley 1994; 1996; Helmi et al. 1999; Chiba & Beers 
2000, hereafter CB; Yanny et al. 2000) and of the Sagittarius dwarf galaxy, which is presently 
being disrupted by the Galactic tidal field (Ibata, Gilmore, h Irwin 1994; Ibata et al. 2000), may 
lend further support to this picture. 

Although halo formation via hierarchical assembly of subgalactic systems, such as dwarf 
galaxies, may continue to the present day (Bland-Hawthorn &; Freeman 2000), a large fraction of 
the stellar halo, especially the inner part where the disk lies, should have been completed prior to 
disk formation, since otherwise the disk component is made significantly thicker than is observed 
due to dynamical heating from infalling masses (Toth k. Ostriker 1992). A clear age gap between 
the (thin) disk and the stellar halo supports that the latter consists of ancient populations (e.g. 
Liu &: Chaboyer 2000). Also, studies of star-forming histories in the disk component indicate that 
the disk has been accumulated at an approximately constant rate over the last several billion 
years (Twarog 1980; Sommer-Larsen k. Yoshii 1990); frequent mergings of dwarf galaxies over the 
Galaxy's lifetime will entirely modify the photometric and spectroscopic properties of the disk. 
Thus, one may well postulate that the inner part of the stellar halo, at say i? < 15 kpc, retains a 
fossil imprint of how it was formed. 

The question then arises, what was the structure of the halo component prior to disk 
formation? Since the bulk of halo stars are found in this inner part, where the disk gravity is 
dominant, the present-day structure of the halo can be greatly afi'ected by later disk formation. 
It is thus necessary to consider the dynamical effect of the disk in inferring the structure of the 
ancient halo from the currently observed halo stars. 

Binney k May (1986, hereafter BM) examined this issue by assuming adiabatic invariance 
(for halo stars) of the action integrals of motion J, and for the distribution function /(J) during 
slow disk formation. They set up test particles distributed in a spheroid, similarly to the stellar 
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distributions of elliptical galaxies, and calculated the dynamical response of the particles to the 
slow increase of disk mass inside the spheroid. They showed that the Galactic halo, before the 
disk was formed, may have had a somewhat flattened density distribution (axial ratio q ~ 0.7), in 
order to produce a current highly flattened halo {q ~ 0.3), which they inferred from the radially 
anisotropic velocity ellipsoid of local metal-poor stars. 

We note here that recent kinematic data for larger samples of local metal-poor stars indicate 
a more moderately flattened halo (g ~ 0.7) for inner radii {R < 15 kpc), whereas the outer halo 
is nearly spherical (Sommer-Larsen & Zhen 1990, hereafter SLZ; CB). This is also supported by 
examination of the spatial distributions of other halo tracers (e.g., Hartwick 1987; Preston et al. 
1991; Kinman, Suntzeff, & Kraft 1994; Yanny et al. 2000). Also, the extent to which the formation 
of the disk component flattened the halo depends on the unknown initial velocity distribution of 
halo stars, whereas the initial conditions set up by BM apply only for one speciflc case. Thus, 
it is yet unexplored what the currently available data for metal-poor stars may tell us about the 
structure of the ancient halo before the disk was formed. 

In this paper we revisit this issue, based on a large sample of halo stars in the solar 
neighborhood, taken from a recently completed catalog of metal-poor stars selected without 
kinematic bias (Beers et al. 2000). It is noted that in a similar vein, Sommer-Larsen (1986), in his 
thesis work, arrived at a conclusion similar to BM's by investigating the distribution of individual 
orbital inclinations for his sample of 143 stars with [Fe/H]< —1.2. In contrast, we seek herein 
to develop a more general and direct method, based on the BM picture, to calculate the global 
density distribution of the halo prior to disk formation. The method is then applied to the more 
accurate and numerous data for metal-poor stars that is presently available. 

This paper is organized as follows. In §2 we describe general properties of orbits in the mass 
model of the Stackel type that we adopt here, as well as the methodology for constructing the 
global density of a given local sample both before and after disk formation. The mass model for 
the Galactic potential that we adopt, and the local sample of metal-poor stars used in the current 
analysis, are also described. In §3 we compute the orbital motions of the sample stars while 
conserving action integrals of motion. We then present the adiabatic change of the derived density 
distributions and kinematic properties of the sample when the disk component is slowly removed 
(essentially working backwards from the present to the past). Finally, in §5, the results are 
summarized, and their implications for the formation and evolution of the Galaxy are discussed. 



2. Method and Sample 

In this work, we adopt the axisymmetric Galactic potential of the Stackel type, for which 
the Hamilton- Jacobi equation separates in spheroidal coordinates (e.g., de Zeeuw 1985; Dejonghe 
& de Zeeuw 1988). Although this type of potential omits the resonant orbits and accompanied 
stochastic orbits that can be revealed in a more general non-separable potential, the fraction 
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of such orbits, which are basically associated with the 1:1 resonance in the radial and vertical 
directions, is thought to be small in the Galactic phase space (May & Binney 1986; BM). Moreover, 
in contrast to a non-separable potential, for which extensive numerical integrations of orbits are 
required, the analytic nature of the Stackel type model has the great advantage of maintaining 
clarity in the analysis. The role of resonant and stochastic orbits revealed from the sample we 
investigate below will be discussed elsewhere (Allen et al. 2000). 



2.1. Stackel Models and Orbital Densities 

In the following, we briefly describe the basic properties of the Stackel models and refer the 
reader to, e.g., de Zeeuw (1985) and Dejonghe &; de Zeeuw (1988), for more details. 

We construct the axisymmetric Galactic potential of the Stackel type, which is defined 
in spheroidal coordinates (A, (j), u) , where corresponds to the azimuthal angle in cylindrical 
coordinates {R, (f), z), and A and u are the roots for r of 

+ r^ = i> (1) 



T + a r + 7 

where a and 7 are constants, giving —7 < u < —a < A. The coordinate surfaces are spheroids 

(A = const.) and hyperboloids of revolution (ly = const.) with the z-axis as the rotation axis, 
where the focal distance A = (7 — a)^/^ fixes the coordinate system. 

The gravitational potential of the Stackel type is then written as 

V,(A..) = -<^ + ''"'W-<- + ^)^M, (2) 

A — u 

where G(r) is an arbitrary function. In this work, G(r) is the sum of Gd{t) from a disk and 
Gh{t) from a massive dark halo. 

The Hamiltonian per unit mass, H, for motion in the potential ■0(A, v) is written as 

where P and Q are the metric coefficients of the spheroidal coordinates, given by 

4(A + a)(A + 7)' ^ 4(i/ + a)(z/ + 7) ' ^' 
and px., p^, and are the conjugate momenta to A, cp, and ly, respectively, 

Px = -P^A = Pvx, pcj, = R^^ = Rvtj>, pu = Q^v = Qvj, . (5) 

The velocities v\, v^j,, and v,y at a point {X,(f),u) are the components of the velocity v along the 
orthogonal axis defined locally by spheroidal coordinates. 
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The "standard" three integrals of motion, I = {E, 12,13), are defined as 

E = -H (6) 
h = f (7) 



2 ^ A-i/ 



(8) 



Another useful set of integrals in the current study are the action integrals, J = (Jx, J^, J^), which 
are defined as 

1 r 2 /■'^2 

Jx = ^fP>^dX = - pxdX (9) 

= ^ ^ P<j>d(t) = Lz (10) 

If 2 fo 

Jv = 2^ f P"^^ ^ J Pi^du , (11) 

where (Ai, A2) and are the turning points of the orbit, defined as the values for which vx = 
and Vi, = 0, respectively, and u = —7 defines the equatorial plane. For the evaluation of Jx, we 
have taken four times the integrals from Ai to A2, to maintain symmetry between Jx and J^ and 
ensure continuity of the actions across transitions from one orbital family to another (de Zeeuw 
1985). 

With a set of three integrals of motion, 1= {E, 12,13), a distribution function /(x, v) can be 
expressed as f{E, 12,13) by application of Jeans' theorem, provided that the system has reached 
dynamical equilibrium. Note that we expect this description to be valid in the inner halo of 
the Galaxy, where the dynamical effect of later accreting materials is small. We employ here 
an ensemble of metal-poor stars which act as tracers moving in the Galactic volume. In such a 
discrete case, consisting of A'" stars, / is written as (Statler 1987), 

N 

f = Y^ Cid{E - Ei)d{l2 - hMh - h,i) , (12) 

i=l 

where Cj is the orbit weighting factor to be determined as described below. 
The density p{x) at any point x is then given as, 

TV 

P{^) = X! (^iPorb{Ei, l2,i, l3,i\ x) , (13) 
1=1 

where pcn-b is the density of a single orbit 

2\/^ 1 

pUE,h.h;^) = -ir viV(A)(A + 7)V-iV(.)(^ + 7)V/2 ' ^''^ 
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with 



iV(r) = G(r)- 



I2 



T + a r + 7 



-E . 



(15) 



With the aid of Jr = Jt-{E, I2, 13), a distribution function can also be expressed as 
f{J\,J,j,,Jiy)- When we define / in action space as 



N 



i=l 



where c| is the orbit weighting factor in action space, the density p is now given as 



N 



(16) 



(17) 



i=l 



where is the density of a single orbit for a given J. We note that p*^^^ cannot generally be 
expressed in an analytic form, as was done in eq.(14) for porhi so we will use porb to calculate the 
density distribution of the stellar halo in this study. 

Comparing eq.(13) with eq.(17) and noting p*^^ = \dl/dlS\porb^ we obtain the relation between 
the amplitudes and q as (Statler 1987) 

^ _ i9(Ja, J(j„ Ju) 



(18) 



2.2. Construction of the Global Density of the Current Stellar Halo 

The formulation given above indicates that, once the orbit weighting factor q is determined, 
we can construct the global density of the halo for a given set of orbits, using eqs.(13)-(14). For 
this purpose we follow the strategies argued by May &: Binncy (1986), as implemented in the 
maximum likelihood approach developed by SLZ, as explained below. 

In the case of a continuous distribution function, /(J), the probability at x of a star with 
actions J, or equivalently, the normalized density distribution of its orbit, is given by (May & 
Binney 1986) 

dJ 



-Por6(x|J) - 



^SPorb ■ (19) 



Then, from Bayes' theorem, the probability that a star found at x has actions in the range 6^3 
centered on J is given by 

/P„,,(X|J)/(J)53J ■ ^ ) 

In the discrete case, consisting of N stars at positions xi,...,xjv, with a distribution function given 
in eq.(16), the probability given above leads to (SLZ) 

p _ '^|Porb,i(^j) _ CiPorb,i{x.j) ^^^-^ 
I^k=lClP*orb,ki'^j) Ilk=lCkPorbA^j) ' 
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Then, Cj is determined by maximizing the probability that the star found at Xj=i is on orbit i = 1, 
the star found at Xj=2 is on orbit i = 2, and so forth. The hkehhood function to be maximized is 



where Porb is estimated from eq.(14). Thus, given a set of (xj, Vj), i = 1,...,N, we can calculate the 
integrals of motion {Ei, l2,i, Is^i) from eqs.(6)-(8), evaluate Cj from eq.(22) with porb from eq.(14), 
and then obtain the global density distribution from eq.(13). 

2.3. Method for Recovering the Ancient Halo before Disk Formation 

We now extend the above method to obtain the global density distribution of the halo prior 

to disk formation, under the assumption that the Galactic disk has been formed slowly compared 
with the dynamical timescale of the system, often taken to be on the order of 10^ yr. This 
assumption is well supported by various lines of evidence (Twarog 1980; Sommer-Larsen & Yoshii 
1990; Rocha-Pinto et al. 2000). In such a case, both the action integrals J [eqs.(9)-(ll)] and the 
distribution function /(J) [eq.(16)] are invariant, whereas I and /(I) change accordingly. 

The method is summarized as follows: (1) At the current epoch, the gravitational potential is 
composed of a disk and a dark halo, G = Gd + Gh, and with this potential, we compute the N 
sets of integrals {Ei, l2^i, I^^i) and {Jx,i, J^^i, Jv,i)-, the orbit weighting factors Cj, and the Jacobians 
\d3/dl\i. We then determine the orbit weighting factors, c|, which are adiabatic invariants. (2) 
At the epoch prior to disk formation, the gravitational potential is supposed to be provided by a 
dark halo alone, G = Gh- With this potential, we search for new N sets of integrals (E'^', /2 j, 13 J 
under the condition that the action integrals {Jx^i, J^^i, Jiy,i) for each star are conserved. We then 
calculate the Jacobians \d3/dV\i and estimate the new orbit weighting factors using eq.(18). (3) 
With {E^, /2 i! -^3 j) and and the dark halo potential, we obtain the global density of the stellar 
halo using eq.(13). 

We note that among the new sets of integrals {E'^, j, I's^i), I2 = h because I2 is expressed as 
h = J^/2 from eq.(7) and eq.(lO). The search of the integrals {E'^, Ig J for a given set of the action 
integrals J requires numerical procedures, except for some specific forms of the potential (see e.g., 
Evans, de Zeeuw, & Lynden-Bell 1990). In our experiment using the sample described below, the 
action integrals J\ and Ji, are well-conserved, within a precision of 0(10^^) in their values. 



We now construct the Galaxy model, consisting of a disk and a dark halo, where both 
components are of the Stackel type. Among the existing Galaxy models of the Stackel type in the 
literature, we adopt the model originally constructed by Dejonghe & de Zeeuw (1988) and later 




(22) 



2.4. The Galaxy Model 
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elaborated by Batsleer & Dejonghe (1994, hereafter BD), which takes a Kuzmin-Kutuzov potential 
for both a highly flattened disk and a nearly spherical halo, 

GgraykM Ggravjl - k)M 

^D[r) = Gh[t) = — — — , {^i) 

+ V -7 Vt + 6 + V-7 + 

where M is the total mass, and k is the ratio of the disk mass to the total mass. The parameter 
h is given so as to let this two-component model remain of the Stackel type. Note that when we 
use ao and cd, instead of a and 7, to define the coordinate system for the disk component, as 
a = —a?j-) and 7 = — c|,, the corresponding quantities for the halo component, an and ch-, must be 
set as = a|) + 6 and = c^j + b. 

To make this mass model resemble the real Galaxy, we set the following conditions: (1) the 
circular velocity around the Galactic center Vrot is nearly constant at ~ 220 km s~^ beyond i? ~ 4 
kpc, (2) the local mass density at the Sun {Rq.Zq) = (8.5,0) kpc is 0.1 ~ 0.2 M© pc"^ (BahcaU 
1984; Bahcall et al. 1992), and (3) the surface mass density at the solar radius is about 70 Mq 
pc~^ within z < 1.1 kpc from the plane (Kuijken &: Gilmore 1991; Bahcall et al. 1992). After 
some experimentation, the following parameter values are adopted as a standard case: cn = 0.052 
kpc, CD/an = 0.02, ch = 17.5 kpc, cj^/a^ = 0.99, M = 10^^ j^^^ ^ q q9_ xhus the disk 

mass IS 9 X 10^° Mq, comprising 9 % of the total mass. Figure 1 shows a rotation curve derived 
from the adopted parameter^, demonstrating that the rotation curve is nearly flat beyond i? ~ 4 
kpc, where Vrot{R&) — 220 km s~'^. The local mass density at the Sun is 0.16 Mq pc~^ and the 
surface mass density sX R = Rq is 68 Mq pc~^ within z < 1.1 kpc from the plane. The axial 
ratios of the potential, q^, obtained from these parameters, are 0.84 at i? = 10 kpc, 0.95 at i? = 20 
kpc, and nearly 1 at larger radii. If we remove the disk component from the potential, is very 
close to 1 at all radii. 



2.5. The Sample of Local Metal-Poor Stars 

To construct the global density of the stellar halo based on the orbits of local metal-poor 
stars, it is important to avoid kinematic bias in the selection of the sample. If the halo stars are 
selected based on high proper motions, for example, the sample will have a bias against stars 
which show similar orbital motions to the Sun, and the global density constructed from such a 
sample will also be biased. It is similarly important to use a large and homogeneously analyzed 
sample, to minimize statistical fluctuations in the derived density distribution, and to avoid, or at 
least minimize, other systematic errors. 



^These parameters are basically the same as those adopted by BD, except that the rotation curve in our model 
is approximately flat to well beyond i? ~ 20 kpc, whereas the BD model shows a falling rotation curve beyond this 
radius. We note that in the BD model, some number of the sample stars we will use below are unbound to the Galaxy 
because of the insufficient mass at outer radii. In our model, all sample stars are bound to the Galaxy, although the 
mass distribution beyond i? ~ 20 kpc is essentially irrelevant to our modeling of the stellar halo inside this radius. 
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Beers et al. (2000) presented a large catalog of metal-poor stars with [Fe/H]< —0.6, selected 
without kinematic bias. A subset of 1214 stars in the catalog contain accurate proper motions 
taken from recently completed proper motion catalogs, including the Hipparcos catalog (ESA 1997), 
in addition to other homogeneously analyzed data with updated stellar positions, newly-derived 
homogeneous distance estimates, revised radial velocities, and refined metal-abundance estimates. 
This is by far the largest sample of metal-poor field stars with available proper motions among 
extant non-kinematically selected samples. Thus, the sample is most advantageous for our current 
study. 

We select, as representatives of the halo population, the stars in the sample within the 
abundance range [Fe/H]< —1.5, which is sufficiently metal-poor to avoid contamination from stars 
with disk- like kinematics (Chiba & Yoshii 1998; CB). Also, to minimize the effects of distance 
errors, we confine ourselves to the stars with measured distances D < 2.5 kpc and with rest-frame 
velocities < 550 km s~^, which is a likely range to bind stars inside the Galaxy. The latter limit 
excludes only two stars. Furthermore, in order to investigate whether there is a finite difference 
between the density distributions of more metal-poor and metal-rich halo populations, as argued 
by Sommer-Larsen (1986), we arbitrary split the sample into two abundance ranges, [Fe/H]< —1.8 
and —1.8 <[Fe/H]< —1.5, as a standard case. The effect of changing the abundance intervals on 
the result will be given in the later section. After applying these cuts, the sample we investigate 
includes = 321 stars for [Fe/H]< -1.8, and N = 182 stars for -1.8 <[Fe/H]< -1.5. Local 
kinematics of the sample in the solar neighborhood are characterized by radially anisotropic 
velocity dispersions, {crB,o-(i,,az) = (153 ± 6, 115 ± 5,97 ± 4) km for [Fe/H]< —1.8 and 
(147 ± 8,114 ± 6, 81 ± 4) km s"^ for -1.8 <[Fc/H]< -1.5, and slow systematic rotation, 
< Vtf) >= 31 ± 6 km s~^ and 37 ± 8 km s~^, for the respective abundance ranges. 

3. Results 

In this section we investigate the adiabatic change of orbital properties and density 
distributions of our sample stars when the disk mass is removed adiabatically from the total 
Galactic potential. 

3.1. Adiabatic Change of the Individual Orbits 

For each star in the potential, after and before disk formation, we compute apo- and 
peri-Galactic distances along the Galactic plane {Rap,k, Rpr,k), and the maximum height away 
from the plane Zmax,ki where k = I for the current potential with the disk and k = 2 when 
the disk is removed adiabatically. In Figures 2a and 2b we show the change of these distances, 
demonstrating that before disk formation, the halo stars orbit at systematically larger distances, 
and similarly, the spatial ranges of the pre-disk orbits are larger than for the current epoch. This 
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is explained by the fact that when a disk is slowly formed in the inner region of the proto-Galactic 

sphere, the gravitational potential becomes more centrally concentrated, and both of the integrals 
{E,!^) are accordingly increased, so that the allowed regions for orbital motions are reduced and 
shifted toward the Galactic center. The characteristic expansion factor in Rap, when the disk is 
removed, is estimated as Rap,2/Rap,i — 1-3. It is also noted that the change of the potential is 
more pronounced in the z-direction because of the disk geometry. As a result, orbital inclinations 
with respect to the plane will be reduced after disk formation (Yoshii & Saio 1979). This is 
actually observed in our calculations, as demonstrated in Figure 2c, where we plot the change of 
inclination angles as defined by ^ = aictaii{zmax / Rap) ■ 

We also compute the orbital eccentricities, defined as e = {rap — rpr) / {vap + rpr) , where rap and 
rpr denote apo- and peri-Galactic distances from the Galactic center, respectively, and plot them 
in Figure 2d. It is apparent that the orbital eccentricities derived here show only a little change, 
especially for e < 0.4 and e > 0.8, even though both rap and rpr change substantially. Thus, the 
orbital eccentricities, defined arbitrarily as above, are approximately adiabatic invariants during 
slow disk formation, as was also demonstrated by Eggen, Lynden-Bell, Sandage (1962) and Yoshii 
& Saio (1979). 

3.2. Global Density Distributions of the Stellar Halo 

Following the method outlined in §2, we now calculate the global density distributions of the 
stellar halo, both at the current epoch and before disk formation. As was done by SLZ and CB, 
we proceed to average the density distributions derived from eq.(13) over grids of finite area in 
the meridional plane of the spheroidal coordinates, (A, v). The grids are defined as Xk = — a, 
k = l, ...,30 and f , = (7 - a) cos2(6';) - 7, 6*/ = (7r/2)(Z/20), I = 0, ...20. The spatial resolution of 
the grids is about 1 kpc. 

In Figure 3a we plot, for the [Fe/H]< —1.8 sample, the radial density distributions along the 
Galactic plane (the averaged density over the area at / = 20) at the current epoch (open circles), 
and when the disk is removed adiabatically (filled circles). As was shown by SLZ and CB, the 
density distribution for R > 8 kpc is well described by a power-law model p oc R^. At the current 
epoch, we find an exponent /3 = —3.4 over 8 < < 30 kpc, in good agreement with the results by 
SLZ and CB. Below R = 8 kpc, the density distributions clearly deviate from a single power-law 
model, a result which is likely caused by incomplete representation of stars with apocentric radii. 
Rap, below R ~ Rq (SLZ; CB) in the local samples we investigate. When the disk is removed 
adiabatically, the density distributions arc made shallower: we obtain an exponent /3 = —3.0 over 
10 < i? < 30 kpc, where the lower radius for this estimate, 10 kpc, is increased from 8 kpc by 
taking into account the characteristic expansion factor of Rap obtained in the previous subsection. 
Thus, as expected, the density distribution of the stellar halo is made more centrally concentrated 
when the disk is slowly formed in the central region of the dark halo. For —1.8 <[Fe/H]< —1.5, 
power-law models with exponents P = —3.3 at the current epoch, and with /3 = —3.0 when the 
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disk is removed, provide excellent fits to the data (Figure 3b). 

In order to obtain a typical error in the estimate of the exponent /?, which arises from the 
combined effects of observational errors in positions and velocities of the sample stars, we have 
constructed ten sets of "pseudo-data" for positions and velocities, where each value is randomly 
selected within its standard observational error with respect to its mean value. From independent 
analysis of these ten reconstructed models, wc find a rms error of 0.14 in the determination of (3. 

Figures 4a and 4b show the equidensity contours of the constructed global density distributions 
in the {R,z) plane for [Fc/H]< —1.8 and —1.8 <[Fc/H]< —1.5, respectively, at the current epoch 
(left panels) and before disk formation (right panels). The lack of stars at small R and large 
z (which gives rise to the ill-formed contour levels in this portion of the diagram) is again a 
consequence of the small probability that stars in the Galaxy which explore such a region are 
represented in the solar neighborhood, as argued in SLZ and CB. This is also seen in Figure 2c, 
where the inclination angles of orbits with respect to the plane, (, are mostly confined to less 
than about 45°, as is the global density inferred from such orbits. Excluding this region, these 
equidensity contours suggest clearly that the current density distributions are flattened at inner 
radii and round at outer radii, as was obtained by CB using a different Galactic potential. In 
contrast, when the disk is removed adiabatically, the density distributions are made rounder, 
especially at inner radii. 

To be more quantitative, we fit elliptical contours to the constructed density maps while 
excluding the region with polar angle 6 < 45°. Specifically, we obtain fits to ellipses of major axis 
a and axial ratio q. The change of our estimate of q' as a function of radius is shown in Figures 5a 
and 5b for [Fe/H]< —1.8 and —1.8 <[Fe/H]< —1.5, respectively, where the error bars correspond 
to the rms errors from the best fits. 

First, at the current epoch (open circles), the axial ratio q remains small at i? < 15 kpc 
and increases with R at larger radii, in good agreement with CB. For [Fe/H]< —1.8, we obtain 
q ~ 0.70 - 0.75 at < 15 kpc and q ~ 0.95 at i? = 20 kpc. For -1.8 <[Fe/H]< -1.5, we obtain 
q ~ 0.50 — 0.60 at i? < 15 kpc and g 2± 1.0 at i? = 18 kpc. The latter subsample appears to show 
a dip in g at i? ~ 20 kpc, which wc believe is a statistical fluctuation due to the limited size of 
the sample (N = 182). It should be noted that, at i? ~ 20 — 25 kpc, the tangential anisotropy of 
the velocity dispersions begins to dominate (Sommcr-Larscn et al. 1997), so that the probability 
that the stars which explore to large z are represented in the solar neighborhood may be small 
for R > 20 kpc. Except for such large radii, it is interesting that the more metal-rich halo 
subsample shows a more flattened density distribution at inner radii than the more metal-poor 
halo subsample. This is consistent with the somewhat smaller vertical velocity dispersion cr^ at 
R = Rq, in the former (81 km s~^) as compared to the latter (97 km s~^). 

Second, when the disk is removed adiabatically (filled circles), the axial ratio is increased to 
q ~ 0.80 for [Fe/H]< -1.8 and q ~ 0.70 for -1.8 <[Fe/H]< -1.5. This result indicates that the 
density distribution of the stellar halo in the inner part of the halo is rounder before disk formation. 
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It should be noted here that this explanation is valid in the global sense: for instance, at i? ~ 20 

kpc in Figure 5a, the axial ratio is actually decreased before disk formation, which is caused by 
the radial expansion of the inner, more flattened region by removing the disk component. The 
result also indicates that the axial ratio q before disk formation, in the assumed spherical potential 
with = 1, is smaller than q^, possibly due to the anisotropic velocity ellipsoid of the stars even 
before disk formation (see van der Marel 1991 for examples of the comparison between q and q^). 
Also, there is a tendency of a more flattened density distribution for —1.8 <[Fc/H]< —1.5 than for 
[Fe/H]< —1.8 even before disk formation, although its significance is too small to be certain. 

In addition to the standard case shown above, we investigate various cases by changing the 
model parameters for the potential or the conditions for the sample selection. Table 1 shows 
the axial ratios of the derived density distributions after (qi) and before (^2) disk formation, 
where, for the sake of straightforward comparison, we have calculated an average axial ratio over 
8 < i? < 13 kpc for qi and 8 < i? < 17 kpc for q2. First, when the disk mass is made smaller 
than the standard case {k = 0.09 0.08), the present-day potential becomes more spherical, 
and so qi is increased. The amount of the change from qi to q2 is the same as or somewhat 
smaller than the standard case, while a small difference in q2 between two abundance ranges still 
remains. Second, a similar change of this result is also seen when we employ the more flattened 
dark halo {cn/dH = 0.99 0.98). In this case, the contribution of the dark halo to the midplane 
potential is increased at a given radius, in a manner that the shape of the potential is made more 
spherical than the standard CclSGj clS would pertain to a less massive disk. Third, we adopt the 
more restrictive distance cut (D < 2.5 kpc — > 1 kpc), to eliminate distant stars for which the 
errors in the estimated velocities are generally larger. It follows from Table 1 that while the result 
for [Fe/H]< —1.8 remains essentially the same, both qi and q2 are systematically decreased for 
— 1.8 <[Fe/H]< —1.5: the difference in q2 between two abundance ranges is increased. Fourth, 
we change the abundance ranges so that the boundary value dividing the two ranges is decreased 
([Fe/H]rfj^ = —1.8 —1.9 or —2.0). Since the more metal-poor stars are included in the more 
metal-rich subsample as [Fe/Hj^i^ is decreased, the difference in the ax;ial ratios turns out to be 
reduced compared to the standard case. 

To summarize, the currently flattened density distribution of the stellar halo in the inner 
region is due both to adiabatic flattening caused by the slowly formed disk and to an initially 
(slightly) flattened density distribution. The tendency that the more metal-rich halo sample 
exhibits a larger flattening even before disk formation may offer a clue for understanding the 
formation process of the stellar halo (see §4). 
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3.3. Adiabatic Change of the Velocity Dispersions 

We now compute the velocity dispersions of the stars at any point x predicted by the current 
model (following Dejonghe k. de Zeeuw 1988 and SLZ), 



N 



^ ^ ^lZciPor6(x)u?,j(x) , r = A,i/ (24) 



^2 1 



N 



1 ^ 



where the velocities {vx, u^, Vj/) are defined by 



PW i=i 



(25) 




^^ = Ml X ^ ^<^ = ±^> = \ J' , (26) 
with 

It = {X + ^)[G{X)-E]-^h (27) 

I3 = {i^ + l)[G{u)-E]-^l2. (28) 

ly + a 

Figure 6a shows the predicted velocity dispersions of the stars with [Fe/H]< —1.8, in the 
Galactic plane {z = 0) at radii of i? = 8.5 and 15.5 kpc. The abscissa denotes the axial ratio of 
the potential, g^, at these radii: = 0.80 at i? = 8.5 kpc and 0.93 at i? = 15.5 kpc when the disk 
is in place, and = 1.00 at both radii before disk formation. Each line shows the change of cJt- 
(r = A, (j), u) after and before disk formation, with endpoints drawn as filled circles for R = 8.5 kpc 
and open circles for R = 15.5 kpc. For instance, at i? = 8.5 kpc, {(T\, a^, a^) = (158, 113, 112) km 
s~^ at the current epoch (left-hand filled circles) and (107,83,64) km s^^ before disk formation 
(right-hand filled circles). We note here that based on the ten pseudo-data sets described in §3.2, 
we find rms errors of these velocity dispersions as (4,3,4) km s"^ at i? = 8.5 kpc and (16,3,4) 
km s~^ at i? = 15.5 kpc. Analogous to panel (a). Figure 6b shows the ratios of these velocity 
dispersions {(J(f)/crx,au/ax) at the same radii. These figures, which arc to be compared with 
Figures 4 and 5 in BM, indicate that the velocity dispersions increased as the disk was formed, as 
a consequence of a more centrally concentrated potential arising from the disk component. Also, 
as panel (b) shows, crj, is boosted more readily than due to the flattening of the potential, so 
that GvIg\ increases. This change of is more prominent at inner radii, where the disk mass 
density is large. On the other hand, the ratio o^jox decreases as the disk was formed. This 
may be explained in the following way. The epicycle theory of orbits (e.g., Binney & Tremaine 
1987) shows that the ratio of velocity dispersions along the plane, a^/ap-,-, is approximately equal 
to K/2fl, where k and Q, denote epicyclic and rotational frequencies, respectively, at a given 
radius. Thus, since k/2CI decreases if the mass distribution responsible the potential is made more 
centrally concentrated, a^/aR (or equivalently a^/ax near the plane) decreases. 
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Although these results are basically consistent with the general trend of the BM model, there 
are notable differences in the values obtained. BM showed that their ag along the polar angle is 
only 38 km at R = Rq before disk formation, suggesting significant anisotropic flattening of 
the initial density distribution {q ~ 0.7). In contrast, the corresponding value using our ai, is 64 
km and the density is more moderately flattened {q 2± 0.8). Thus, the velocity dispersions 
of the halo component before the disk was formed are characterized by an anisotropic velocity 
ellipsoid, but its degree of anisotropy is more moderate than previously obtained. 

The more metal-rich abundance range, —1.8 <[Fe/H]< —1.5, gives basically the same 
result as the above, except for the ratios of the velocity dispersions before disk formation - 
{a^/(Tx,(Tu/(Tx) at R = 8.5 kpc is (0.59,0.44), for -1.8 <[Fe/H]< -1.5, but takes values (0.77,0.60) 
for [Fc/H]< —1.8. Thus, the former, more metal-rich range, shows a more anisotropic velocity 
ellipsoid than the latter. This is consistent with the more flattened density distribution for the 
former, even before disk formation, as was shown in the previous subsection. 

4. Discussion and Conclusion 

The global structure of the present-day stellar halo is characterized by an inner, highly 
flattened part, as revealed at < 15 kpc, and an outer, nearly spherical part (SLZ; CB). This 
two-component picture for the present-day stellar halo provides a reasonable explanation why 
faint-star-count studies have generally yielded an approximately spherical halo, whereas the local 
anisotropic velocities of the halo stars suggest a highly flattened system (Freeman 1987). The 
issue relevant here is what physical mechanism in the early stage of the Galaxy gives rise to the 
inner, highly flattened halo, where the bulk of halo stars are found, and where the effects of later 
satellite accretion may be diminished. 

One of the possible reasons for the two-component nature of the present-day halo is the 
slow formation of the disk within an initially spherical stellar halo (BM). In this paper, we have 
quantified the effect of later disk formation on the halo flattening, based on methods assuming 
adiabatic invariance of the motion of halo stars, and its application to a large sample of stars in the 
solar neighborhood. We have found that, even before disk formation, the inner part of the stellar 
halo exhibited a finite flattening, although it is more moderate than presently observed. The axial 
ratios of the density profiles within the almost spherical potential are q ~ 0.80 for [Fe/H]< —1.8 
and q ~ 0.70 for —1.8 <[Fe/H]< —1.5. Also, the initial velocity dispersions are characterized by 
an anisotropic velocity ellipsoid, as {(Tci,/a\,ay/ax) = (0.77,0.60) for [Fe/H]< —1.8 and (0.59,0.44) 
for —1.8 <[Fc/H]< —1.5 at i? = Rq. Therefore, the inner part of the stellar halo was flattened by 
velocity anisotropy, and the more metal-rich population likely exhibited a more flattened density 
distribution. 

Through a comparison of the inclination angles of orbits for two abundance ranges, 
[Fe/H]< —1.5 and —1.5 <[Fe/H]< —1.2, Sommer-Larsen (1986) also obtained a more flattened 
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initial distribution for more metal-rich populations. We note that his latter subsample shows a 
rapid systematic rotation of 147 km s~^, possibly contaminated by the stars belonging to the 
metal- weak thick disk (Freeman 1987). In contrast, our metal-rich halo subsample is selected from 
the more restrictive range —1.8 <[Fe/H]< —1.5, where the effect of disk-like kinematics is minimal 
(Chiba &; Yoshii 1998; CB) - the flattening of this subsample is caused by velocity anisotropy, not 
systematic rotation. 

The results presented here may suggest that the ancient halo, at least in its inner part, may 
have undergone a somewhat ordered contraction. This contraction could have involved dissipation 
due to baryonic gas - radiative cooling of this gas was most efficient in the innermost regions with 
high density. The resultant contraction of this gas, as a whole, may have proceeded mainly along 
the axis of rotation, because of the absence of the angular-momentum barrier in this direction. As 
the chemical enrichment proceeded along with the progress of the collapse, more metal-enriched 
stars would have "seen" a more flattened density distribution. On the other hand, the outer part 
of the halo may have been more susceptible to later infall of satellite galaxies, so that its density, 
kinematics, and mean age are different from those of the inner halo (Norris 1994; Carney et al. 
1996; Sommer-Larsen et al. 1997). 

Alternatively, one might argue that the two-component nature of the present-day stellar halo 
is entirely a consequence of satellite accretion after disk formation (Freeman 1987). According 
to Quinn & Goodman (1986) (see also Quinn, Hernquist, & FuUagar 1993), prograde satellite 
orbits that are initially inclined at less than about 60° to the disk are dragged down quickly 
toward the plane by the effects of the dynamical friction of the disk. Walker, Mihos, &: Hernquist 
(1996) further explore the effects of mergers of small satellites with large disk galaxies such as the 
Milky Way. The debris from these merging satellites, in combination with disrupted disk stars, 
would be expected to form a flattened system. Furthermore, if more massive satellites were more 
metal-rich (as they appear to be at present, see Mateo 1998), their orbits would fall farther toward 
the Galactic center so that their debris would form a more flattened, more metal-rich system 
(Freeman 1987). However, there exists no clear evidence for the predicted dynamical heating of 
the thin-disk component - its very thin geometry (Toth & Ostriker 1992), and the nearly constant 
velocity dispersion of thin-disk stars over the last 10 Gyrs (Quillen & Garnett 2000), suggest that 
the thin disk has sustained little significant damage since its formation. In this regard, one might 
argue that the metal-weak thick disk is evidence for early dynamical heating of a pre-existing, 
metal-deficient, thin disk. However, it is then difficult to explain its absence in the abundance 
range considered in this paper (Chiba & Yoshii 1998; CB). 

It is worthwhile to remark that the hypothesis of the dissipative formation of the inner 
flattened halo, as well as the later accretion of satellites onto the outer halo, is a natural 
consequence of the CDM hierarchical clustering model (Bekki & Chiba 2000). This model 
postulates that a protogalactic system initially contains numerous subgalactic clumps, comprised 

of a mixture of gas and dark matter, and that the merging of these clumps led to a smaller number 
of more massive clumps. In the simulations of Bekki & Chiba (2000), these larger clumps move 
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gradually toward the central region of the system, due to both dynamical friction and dissipative 

merging with smaller clumps. Finally, the last merging event occurs between the two most massive 
clumps, and the metal-poor stars which have been formed inside the clumps are disrupted and 
spread over the inner part of the halo. The aftermath is characterized by a flattened density 
distribution. Some fraction of the disrupted gas from the clumps may settle into the central region 
of the system, and produce a more enriched, more flattened density distribution. Some of the 
initially small density fluctuations in the outer region would have gained systematically higher 
angular momentum from their surroundings, and then slowly fallen into the system after most 
parts of the system were formed. This may correspond to the process of late satellite accretion, 
contributing primarily to the outer part of the halo. Thus, the reported initial state of the stellar 
halo can be explained, at least qualitatively, in the context of hierarchical clustering scenario. 

An alternative approach for elucidation of the dissipative nature of halo formation is to 
examine the results of recent high-resolution N-body simulations of structure formation based on 
the CDM theory (e.g., Ghigna et al. 1998; Moore et al. 1999; Klypin et al. 1999). Such simulations 
provide the orbital properties of dark matter particles inside virialized dark halos. If the stellar 
halo component in the Galaxy is formed similarly through dissipationless hierarchical assembly, 
the orbits of halo stars prior to disk formation, as derived in the current paper, may follow those 
of dark matter particles. For this purpose, we take the Ghigna et al. (1998) simulation of the 
formation of a cluster, as this is currently the only published one that presents the detailed orbital 
distribution of the simulated particles. They showed that the orbital distribution of the halo 
particles is close to isotropic - circular orbits are rare and radial orbits are common. The average 
ratio of pericentric and apocentric distances, rp^/rap, is equal to or less than 0.20, without showing 
a large variation as a function of the distance from the cluster: the median ratio is approximately 
0.17. On the other hand, as is deduced from Figure 2, the orbits of the halo stars we have derived 
here before the disk was in place are more circular than the simulated dark halo particles, and the 
velocity field is anisotropic: the average value of rpr/vap is 0.29. Thus, we require some additional 
process, possibly dissipative interaction among protogalactic clumps to circularize their orbits, to 
explain the characteristic orbital distribution of the halo stars in the early Galaxy. 

More quantitative conclusions must await more elaborate modeling of the formation of the 
Galaxy over a large mimbcr of possible model parameters. Also, it is necessary to assemble and 
analyze the data of more remote stars, especially those presently found inside the solar radius, 
where our modeling of the halo is incomplete. In this regard, the next generation of astrometric 
satellites, such as FAME and GAIA, will provide highly precise parallaxes and proper motions 
for numerous stars, so that both three dimensional positions and velocities will be available over 
a large fraction of the halo. Also, with these astrometric satellites, we will be able to determine 
the exact mass distributions of the disk and dark halo components, and thus obtain definite 
information on the early Galaxy, using the technique outlined here. Furthermore, in addition 
to the Milky Way, direct identification of halo populations in external disk galaxies may prove 
promising as a way to clarify the global structures of stellar halos and their association with 
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disks and bulges (e.g., Morrison 1999). Such studies should be eagerly pursued with 10m class 
telescopes. 

We are grateful to the anonymous referee for constructive comments on the paper. This work 
has been supported in part by Grants-in-Aid for Scientific Research (09640328) from the Ministry 
of Education, Science, Sports and Culture of Japan. 
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Fig. 1. — Rotation curve for the adopted Galactic potential (solid line). Dotted and dashed lines 
denote the contributions from the disk and dark halo, respectively. 
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Fig. 2. — Orbital parameters of the sample stars before (ordinate, subscript 2) and after (abscissa, 
subscript 1) disk formation, for (a) apogalactic distance along the plane Rap, (b) maximum height 
away from the plane Zmax, (c) inclination angle with respect to the plane ^ = SLrctSLn{zmax / Rap) , 
and (d) orbital eccentricity e. 
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Fig. 3. — Density distributions of the reconstructed halo in the Galactic plane before (filled circles) 
and after (open circles) disk formation, for (a) [Fe/H] < -1.8 and (b) -1.8 < [Fe/H]< -1.5. All 
plots have been shifted arbitrarily along the vertical axis for clarity. The dashed and dotted 
lines denote the best-fit power-law model with labeled exponents, before and after disk formation, 
respectively. 
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Fig. 4. — Equidensity contours for the reconstructed halo in the (R, z) plane, for (a) [Fe/H]< —1.8 
and (b) —1.8 <[Fe/H]< —1.5. Left and right panels correspond to the density distributions after 
and before disk formation, respectively. 



-25- 




R (kpc) 

Fig. 5. — Axial ratios for the density distributions of the reconstructed halo before (filled circles) 
and after (open circles) disk formation, for (a) [Fe/H] < —1.8 and (b) —1.8 <[Fe/H]< —1.5. 
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Fig. 6. — (a) Velocity dispersions of the stars with [Fe/H]< —1.8, {a\,a^,a,^), and (b) their ratios, 
(c^/cajCiz/ca)) at two radii R = 8.5 kpc and 15.5 kpc in the Galactic plane {z = 0). The abscissa 
denotes the axial ratio of the potential, q^, at both radii: q^p = 0.80 at R = 8.5 kpc and 0.93 at 
R = 15.5 kpc when the disk is in place, and q^ = 1.00 at both radii before disk formation. Left and 
right hand filled circles show or ar/cr\ after and before disk formation, respectively, for R = 8.5 
kpc, whereas open circles for R = 15.5 kpc. 



Table 1. Axial Ratios of the Stellar Halo for Various Cases 



Model 




[Fe/H] 


<[Fc/H]rf„ 


[Fc/H]<j„ 


<[Fe/H]< -1.5 


gi" 








Standard Case** 




0.72 


0.80 


0.53 


0.72 


Parameter Change 


k = 0.09 ^ 0.08 




0.80 


0.88 


0.55 


0.68 


cn/aH = 0.99 ^ 0.98 




0.74 


0.86 


0.59 


0.80 


D < 2.5 kpc 1 kpc 




0.78 


0.84 


0.44 


0.56 


[Fc/H]d„ = -1.8 ^ - 


1.9 


0.72 


0.80 


0.57 


0.73 


[Fc/H]d„ = -1.8 ^ - 


2.0 


0.65 


0.76 


0.66 


0.83 



'^Our standard choice of the parameters includes D < 2.5 kpc, [Fe/Hj^ju = 
1.8, k = 0.09 (the disk mass is 9 x lO^" Mq), and cn/aH = 0.99. 
''Average axial ratio over 8 < i? < 13 kpc in the present-day disk+halo system. 
■^Average axial ratio over 8 < R < 17 kpc before disk formation. 



